Heatmaps

Ubiquitous

experimental_metadata <- read.delim("data/ubiquitous/metadata_rem.txt", sep=",", header=TRUE, stringsAsFactors=FALSE)
annotation = data.frame(Condition = rep(c("Control","OKSM", "Senolytic", "Senolyic_OKSM"),
                                            c(3, 3, 2, 3)))

row.names(annotation) = experimental_metadata$sample_id 
anno_colours = list(Condition = c(Control = "#BF3100", OKSM = "#E9B44C", Senolytic = "#1B998B", Senolyic_OKSM = "#5D576B"))
pheatmap(ubiq_pam_clust[,1:(ncol(ubiq_pam_clust)-1)],
         color = colorRampPalette(rev(brewer.pal(n = 7, name = "RdBu")))(100),
         fontsize_row = 5.5,
         annotation_col = annotation,
         annotation_colors = anno_colours,
         cluster_rows = FALSE,
         cluster_cols = FALSE)

ubi_pam_clust_df %>% 
  dplyr::select(gene_name, Cluster) %>%
  datatable(options = list(scrollX = TRUE), class = "white-space: nowrap")
ubi_c1 <- ubi_pam_clust_df[ubi_pam_clust_df$Cluster == 1, ] %>%
  dplyr::select(-Cluster)

ubi_c2 <- ubi_pam_clust_df[ubi_pam_clust_df$Cluster == 2, ] %>%
  dplyr::select(-Cluster)

ubi_c3 <- ubi_pam_clust_df[ubi_pam_clust_df$Cluster == 3, ] %>%
  dplyr::select(-Cluster)

ubi_c4 <- ubi_pam_clust_df[ubi_pam_clust_df$Cluster == 4, ] %>%
  dplyr::select(-Cluster)

ubi_c5 <- ubi_pam_clust_df[ubi_pam_clust_df$Cluster == 5, ] %>%
  dplyr::select(-Cluster)

ubi_c6 <- ubi_pam_clust_df[ubi_pam_clust_df$Cluster == 6, ] %>%
  dplyr::select(-Cluster)

ubi_c7 <- ubi_pam_clust_df[ubi_pam_clust_df$Cluster == 7, ] %>%
  dplyr::select(-Cluster)

data.frame(Cluster = c("Cluster 1", "Cluster 2", "Cluster 3", "Cluster 4", "Cluster 5", 
                       "Cluster 6", "Cluster 7", "Total"),
                        Number_of_genes = c(nrow(ubi_c1), nrow(ubi_c2), nrow(ubi_c3), nrow(ubi_c4),nrow(ubi_c5), 
                                            nrow(ubi_c6), nrow(ubi_c7),
                                            sum(c(nrow(ubi_c1),nrow(ubi_c2),nrow(ubi_c3),nrow(ubi_c4),nrow(ubi_c5),
                                                  nrow(ubi_c6),nrow(ubi_c7)))))
##     Cluster Number_of_genes
## 1 Cluster 1             300
## 2 Cluster 2             201
## 3 Cluster 3             157
## 4 Cluster 4             173
## 5 Cluster 5             147
## 6 Cluster 6             204
## 7 Cluster 7             100
## 8     Total            1282

ISC

pheatmap(isc_pam_clust[,1:(ncol(isc_pam_clust)-1)],
         color = colorRampPalette(rev(brewer.pal(n = 7, name = "RdBu")))(100),
         fontsize_row = 5.5,
         annotation_col = annotation,
         annotation_colors = anno_colours,
         cluster_rows = FALSE,
         cluster_cols = FALSE)

isc_pam_clust_df <- isc_pam_clust %>%
  as.data.frame() %>%
  mutate(gene_name = ensembl.genes[rownames(.),]$external_gene_name)
datatable(isc_pam_clust_df, options = list(scrollX = TRUE), class = "white-space: nowrap")
isc_c1 <- isc_pam_clust_df[isc_pam_clust_df$Cluster == 1, ] %>%
  dplyr::select(-Cluster)

isc_c2 <- isc_pam_clust_df[isc_pam_clust_df$Cluster == 2, ] %>%
  dplyr::select(-Cluster)

isc_c3 <- isc_pam_clust_df[isc_pam_clust_df$Cluster == 3, ] %>%
  dplyr::select(-Cluster)

isc_c4 <- isc_pam_clust_df[isc_pam_clust_df$Cluster == 4, ] %>%
  dplyr::select(-Cluster)

isc_c5 <- isc_pam_clust_df[isc_pam_clust_df$Cluster == 5, ] %>%
  dplyr::select(-Cluster)

isc_c6 <- isc_pam_clust_df[isc_pam_clust_df$Cluster == 6, ] %>%
  dplyr::select(-Cluster)

isc_c7 <- isc_pam_clust_df[isc_pam_clust_df$Cluster == 7, ] %>%
  dplyr::select(-Cluster)

data.frame(Cluster = c("Cluster 1", "Cluster 2", "Cluster 3", "Cluster 4", 
                       "Cluster 5", "Cluster 6", "Cluster 7", "Total"),
                        Number_of_genes = c(nrow(isc_c1), nrow(isc_c2), nrow(isc_c3), nrow(isc_c4),
                                            nrow(isc_c5), nrow(isc_c6), nrow(isc_c7),
                                            sum(c(nrow(isc_c1),nrow(isc_c2),nrow(isc_c3),nrow(isc_c4),
                                                  nrow(isc_c5),nrow(isc_c6),nrow(isc_c7)))))
##     Cluster Number_of_genes
## 1 Cluster 1             773
## 2 Cluster 2             397
## 3 Cluster 3             128
## 4 Cluster 4             566
## 5 Cluster 5             613
## 6 Cluster 6             248
## 7 Cluster 7             242
## 8     Total            2967

Gene Expression Boxplots

Ubiquitous

Cluster 1

cluster_boxplot(ubi_c1, "Cluster 1", "TdTom", "Sen_OKSM", "OKSM")

Cluster 2

cluster_boxplot(ubi_c2, "Cluster 2", "TdTom", "Sen_OKSM", "OKSM")

Cluster 3

cluster_boxplot(ubi_c3, "Cluster 3", "TdTom", "Sen_OKSM", "OKSM")

Cluster 4

cluster_boxplot(ubi_c4, "Cluster 4", "TdTom", "Sen_OKSM", "OKSM")

Cluster 5

cluster_boxplot(ubi_c5, "Cluster 5", "TdTom", "Sen_OKSM", "OKSM")

Cluster 6

cluster_boxplot(ubi_c6, "Cluster 6", "TdTom", "Sen_OKSM", "OKSM")

Cluster 7

cluster_boxplot(ubi_c7, "Cluster 7", "TdTom", "Sen_OKSM", "OKSM")

isc

Cluster 1

cluster_boxplot(isc_c1, "Cluster 1", "Control", "SO", "O")

Cluster 2

cluster_boxplot(isc_c2, "Cluster 2", "Control", "SO", "O")

Cluster 3

cluster_boxplot(isc_c3, "Cluster 3", "Control", "SO", "O")

Cluster 4

cluster_boxplot(isc_c4, "Cluster 4", "Control", "SO", "O")

Cluster 5

cluster_boxplot(isc_c5, "Cluster 5", "Control", "SO", "O")

Cluster 6

cluster_boxplot(isc_c6, "Cluster 6", "Control", "SO", "O")

Cluster 7

cluster_boxplot(isc_c7, "Cluster 7", "Control", "SO", "O")

Gene overlaps

Overall

olaps <- readPNG("output/comparison/overlaps.png")
p <- rasterGrob(olaps, interpolate = TRUE)
ggarrange(p)

What are these 686 genes that overlap?

overlapping_genes <- ubi_pam_clust_df[rownames(ubi_pam_clust_df) %in% rownames(isc_pam_clust_df),]
overlapping_genes %>%
  dplyr::select(gene_name, Cluster) %>%
  rename("Gene Name" = gene_name) %>%
  datatable(options = list(scrollX = TRUE))

Where are they located in the ubiquitous dataset?

overlapping_genes %>%
  dplyr::select(gene_name, Cluster) %>%
  rename("Gene Name" = gene_name) %>%
  group_by(Cluster) %>%
  tally() %>%
  datatable(options = list(scrollX = TRUE), caption = "Where the 686 genes are located in the ubiquitous dataset")

What are these overlapping genes enriched for?

GO

eresList$GO_Biological_Process_2018 %>%
  plot_enrichr("GO_Biological_Process_2018")

datatable(eresList[[1]], options = list(scrollX = TRUE), class = "white-space: nowrap")
terms_to_highlight_go <- eresList$GO_Biological_Process_2018 %>%
  mutate(Annotated = as.numeric(str_extract(as.character(Overlap), "\\d+$")),
         Significant = as.numeric(str_extract(as.character(Overlap), "^\\d+")),
         Ratio = Significant/Annotated) %>%
  dplyr::filter(Term %in% c("DNA replication (GO:0006260)", "nervous system development (GO:0007399)", "mitotic cytokinesis (GO:0000281)", "proteolysis (GO:0006508)", "regulation of cell proliferation (GO:0042127)"))
terms_to_highlight_go %>%
  mutate(Term = gsub("\\([^()]*\\)", "", Term),
         Term = str_to_title(Term)) %>%
  dplyr::arrange(Ratio) %>%
  mutate(Term = factor(Term, levels = Term[rev(order(Ratio))])) %>%
  ggplot(aes(x = Ratio, y = Term)) +
  geom_point(aes(colour = Adjusted.P.value, size = Ratio)) + 
  ylab(NULL) + 
  xlab("Gene Ratio") +
  scale_color_continuous(low="red", high="blue", guide=guide_colorbar(reverse=FALSE)) + 
  scale_y_discrete(labels = function(x) str_wrap(x, width = 60),
                   limits=rev) +
  #labs(fill = "Adjusted p-value") +
  guides(colour = guide_colorbar(title = "Adjusted p-value", reverse = TRUE),
         size = guide_legend(title = "Gene Ratio")) +
  theme(axis.text.x = element_text(size=8)) +
  theme_bw() +
  ggtitle("KEGG Enrichment Analysis")

KEGG

eresList$KEGG_2019 %>%
  plot_enrichr("KEGG_2019")

datatable(eresList[[2]], options = list(scrollX = TRUE), class = "white-space: nowrap")
terms_to_highlight_kegg <- eresList$KEGG_2019 %>%
  mutate(Annotated = as.numeric(str_extract(as.character(Overlap), "\\d+$")),
         Significant = as.numeric(str_extract(as.character(Overlap), "^\\d+")),
         Ratio = Significant/Annotated) %>%
  dplyr::filter(Term %in% c("Toll and Imd signaling pathway", "Apoptosis", "DNA replication", "Hippo signaling pathway", "Lysosome"))
terms_to_highlight_kegg %>%
  mutate(Term = gsub("\\([^()]*\\)", "", Term),
         Term = str_to_title(Term)) %>%
  dplyr::arrange(Ratio) %>%
  mutate(Term = factor(Term, levels = Term[rev(order(Ratio))])) %>%
  ggplot(aes(x = Ratio, y = Term)) +
  geom_point(aes(colour = Adjusted.P.value, size = Ratio)) + 
  ylab(NULL) + 
  xlab("Gene Ratio") +
  scale_color_continuous(low="red", high="blue", guide=guide_colorbar(reverse=FALSE)) + 
  scale_y_discrete(labels = function(x) str_wrap(x, width = 60),
                   limits=rev) +
  #labs(fill = "Adjusted p-value") +
  guides(colour = guide_colorbar(title = "Adjusted p-value", reverse = TRUE),
         size = guide_legend(title = "Gene Ratio")) +
  theme(axis.text.x = element_text(size=8)) +
  theme_bw() +
  ggtitle("KEGG Enrichment Analysis")

Where in the ISC dataset are these 686 genes?

isc_pam_clust_df[rownames(isc_pam_clust_df) %in% rownames(ubi_pam_clust_df),] %>%
  dplyr::select(gene_name, Cluster) %>%
  group_by(Cluster) %>%
  tally() %>%
  datatable(options = list(scrollX = TRUE), caption = "Where the 686 genes are located in the ISC dataset")

I’ve manipulated the heatmaps such that Clusters 1 - 5 in the ubiquitous data is the same as Clusters 1 - 7 in the ISC dataset for easy comparison

C2 ubiquitous dataset in C1 and/or C2 of ISC dataset

Percentage of genes in C2 of ubiquitous dataset in c1 & c2 of ISC dataset

denom <- c(rownames(isc_c1),rownames(isc_c2))
sum(rownames(ubi_c2) %in% denom)/length(rownames(ubi_c2))*100
## [1] 34.32836

Where else can we find c2 ubiquitous genes in ISC dataset?

tab <- isc_pam_clust_df[rownames(ubi_c2),] %>%
  drop_na %>%
  group_by(Cluster) %>%
  tally()
tab
## # A tibble: 7 × 2
##   Cluster     n
##     <dbl> <int>
## 1       1    39
## 2       2    30
## 3       3     2
## 4       4     3
## 5       5     6
## 6       6    12
## 7       7     4

96 genes in Cluster 2 of the ubiquitous dataset is accounted for in the ISC dataset. This implies that 204 genes in C2 of the ubiquitous dataset is not expressed or significantly differentially expressed in the ISC dataset.

C3 ubiquitous dataset in C3 of ISC dataset

Percentage of genes in c3 of ubiquitous dataset in c3 ISC dataset

sum(rownames(ubi_c3) %in% rownames(isc_c3))/length(rownames(ubi_c3))*100
## [1] 11.46497

Where else can we find c3 ubiquitous genes in ISC dataset?

tab <- isc_pam_clust_df[rownames(ubi_c3),] %>%
  drop_na %>%
  group_by(Cluster) %>%
  tally()
tab
## # A tibble: 7 × 2
##   Cluster     n
##     <dbl> <int>
## 1       1    22
## 2       2    28
## 3       3    18
## 4       4     6
## 5       5     5
## 6       6     5
## 7       7     6

90 genes in Cluster 3 of the ubiquitous dataset is accounted for in the ISC dataset. This implies that 67 genes in C3 of the ubiquitous dataset is not expressed or significantly differentially expressed in the ISC dataset.

C4 ubiquitous dataset in C4 & C5 of ISC dataset

Percentage of genes in c4 of ubiquitous dataset in c6 ISC dataset

denom <- c(rownames(isc_c4), rownames(isc_c5))
sum(rownames(ubi_c4) %in% denom)/length(rownames(ubi_c4))*100
## [1] 28.90173

Where else can we find c4 ubiquitous genes in ISC dataset?

tab <- isc_pam_clust_df[rownames(ubi_c4),] %>%
  drop_na %>%
  group_by(Cluster) %>%
  tally()
tab
## # A tibble: 7 × 2
##   Cluster     n
##     <dbl> <int>
## 1       1     8
## 2       2     4
## 3       3     2
## 4       4    31
## 5       5    19
## 6       6     9
## 7       7    14

87 genes in Cluster 4 of the ubiquitous dataset is accounted for in the ISC dataset. This implies that 86 genes in C4 of the ubiquitous dataset is not expressed or significantly differentially expressed in the ISC dataset.

C7 ubiquitous dataset in C6 of ISC dataset

Percentage of genes in c4 of ubiquitous dataset in c6 ISC dataset

sum(rownames(ubi_c7) %in% rownames(isc_c6))/length(rownames(ubi_c7))*100
## [1] 18

Where else can we find c7 ubiquitous genes in ISC dataset?

tab <- isc_pam_clust_df[rownames(ubi_c7),] %>%
  drop_na %>%
  group_by(Cluster) %>%
  tally()
tab
## # A tibble: 6 × 2
##   Cluster     n
##     <dbl> <int>
## 1       1    19
## 2       2     7
## 3       3     1
## 4       4     3
## 5       5    10
## 6       6    18

Summary Enrichments from custom script

Ubiquitous

filenames = list.files("output/ubiquitous", pattern = "_go.rds", full.names = TRUE)
all_go = lapply(filenames, function(x){readRDS(x)})
names(all_go) <- c("c1", "c2", "c3", "c4", "c5", "c6", "c7")
ubi_go <- all_go %>%
  bind_rows(.id = "Cluster") %>%
  group_by(Cluster) %>%
  slice_min(order_by = Adjusted.P.value, n = 10)
ubi_go %>%
  mutate(Term = gsub("\\([^()]*\\)", "", Term),
         Term = str_to_title(Term)) %>%
  mutate(Cluster = factor(Cluster, levels = c(unique(ubi_go$Cluster)))) %>%
  group_by(Cluster) %>%
  mutate(Term = factor(Term, levels = Term)) %>%
  ggplot(aes(x = Cluster, y = Term)) +
  geom_point(aes(colour = Adjusted.P.value, size = Ratio)) + 
  ylab(NULL) + 
  xlab("Cluster") +
  scale_color_continuous(low="red", high="blue", guide=guide_colorbar(reverse=FALSE)) + 
  scale_y_discrete(labels = function(x) str_wrap(x, width = 60),
                   limits=rev) +
  #labs(fill = "Adjusted p-value") +
  guides(colour = guide_colorbar(title = "Adjusted p-value", reverse = TRUE),
         size = guide_legend(title = "Gene Ratio")) +
  theme(axis.text.x = element_text(size=8)) +
  theme_bw() +
  ggtitle("GO Enrichment Analysis")

filenames = list.files("output/ubiquitous", pattern = "_kegg.rds", full.names = TRUE)
all_kegg = lapply(filenames, function(x){readRDS(x)})
names(all_kegg) <- c("c1", "c2", "c3", "c4", "c5", "c6", "c7")
ubi_kegg <- all_kegg %>%
  bind_rows(.id = "Cluster") %>%
  group_by(Cluster) %>%
  slice_min(order_by = Adjusted.P.value, n = 10)
ubi_kegg %>%
  mutate(Term = gsub("\\([^()]*\\)", "", Term),
         Term = str_to_title(Term)) %>%
  mutate(Cluster = factor(Cluster, levels = c(unique(ubi_kegg$Cluster)))) %>%
  group_by(Cluster) %>%
  mutate(Term = factor(Term, levels = Term)) %>%
  ggplot(aes(x = Cluster, y = Term)) +
  geom_point(aes(colour = Adjusted.P.value, size = Ratio)) + 
  ylab(NULL) + 
  xlab("Cluster") +
  scale_color_continuous(low="red", high="blue", guide=guide_colorbar(reverse=FALSE)) + 
  scale_y_discrete(labels = function(x) str_wrap(x, width = 60),
                   limits=rev) +
  #labs(fill = "Adjusted p-value") +
  guides(colour = guide_colorbar(title = "Adjusted p-value", reverse = TRUE),
         size = guide_legend(title = "Gene Ratio")) +
  theme(axis.text.x = element_text(size=8)) +
  theme_bw() +
  ggtitle("KEGG Enrichment Analysis")

ISC

filenames = list.files("output/isc", pattern = "_go.rds", full.names = TRUE)
all_go = lapply(filenames, function(x){readRDS(x)})
names(all_go) <- c("c1", "c2", "c3", "c4", "c5", "c6", "c7")
isc_go <- all_go %>%
  bind_rows(.id = "Cluster") %>%
  group_by(Cluster) %>%
  slice_min(order_by = Adjusted.P.value, n = 10)
isc_go %>%
  mutate(Term = gsub("\\([^()]*\\)", "", Term),
         Term = str_to_title(Term)) %>%
  mutate(Cluster = factor(Cluster, levels = c(unique(isc_go$Cluster)))) %>%
  group_by(Cluster) %>%
  mutate(Term = factor(Term, levels = Term)) %>%
  ggplot(aes(x = Cluster, y = Term)) +
  geom_point(aes(colour = Adjusted.P.value, size = Ratio)) + 
  ylab(NULL) + 
  xlab("Cluster") +
  scale_color_continuous(low="red", high="blue", guide=guide_colorbar(reverse=FALSE)) + 
  scale_y_discrete(labels = function(x) str_wrap(x, width = 60),
                   limits=rev) +
  #labs(fill = "Adjusted p-value") +
  guides(colour = guide_colorbar(title = "Adjusted p-value", reverse = TRUE),
         size = guide_legend(title = "Gene Ratio")) +
  theme(axis.text.x = element_text(size=8)) +
  theme_bw() +
  ggtitle("GO Enrichment Analysis")

filenames = list.files("output/isc", pattern = "_kegg.rds", full.names = TRUE)
all_kegg = lapply(filenames, function(x){readRDS(x)})
names(all_kegg) <- c("c1", "c2", "c3", "c4", "c5", "c6", "c7")
isc_kegg <- all_kegg %>%
  bind_rows(.id = "Cluster") %>%
  group_by(Cluster) %>%
  slice_min(order_by = Adjusted.P.value, n = 10)
isc_kegg %>%
  mutate(Term = gsub("\\([^()]*\\)", "", Term),
         Term = str_to_title(Term)) %>%
  mutate(Cluster = factor(Cluster, levels = c(unique(isc_kegg$Cluster)))) %>%
  group_by(Cluster) %>%
  mutate(Term = factor(Term, levels = Term)) %>%
  ggplot(aes(x = Cluster, y = Term)) +
  geom_point(aes(colour = Adjusted.P.value, size = Ratio)) + 
  ylab(NULL) + 
  xlab("Cluster") +
  scale_color_continuous(low="red", high="blue", guide=guide_colorbar(reverse=FALSE)) + 
  scale_y_discrete(labels = function(x) str_wrap(x, width = 60),
                   limits=rev) +
  #labs(fill = "Adjusted p-value") +
  guides(colour = guide_colorbar(title = "Adjusted p-value", reverse = TRUE),
         size = guide_legend(title = "Gene Ratio")) +
  theme(axis.text.x = element_text(size=8)) +
  theme_bw() +
  ggtitle("KEGG Enrichment Analysis")

Summary Enrichments from compareCluster

Ubiquitous

genes = lapply(list(c1=ubi_c1, c2=ubi_c2, c3=ubi_c3, c4=ubi_c4, c5=ubi_c5, c6=ubi_c6, c7=ubi_c7), function(i){ensembl.genes[rownames(i)]$entrezgene_id})
bg = ensembl.genes[rownames(ddssva)]$entrezgene_id
compPathway <- compareCluster(geneCluster   = genes,
                              fun           = "enrichGO",
                              keyType       = "ENTREZID",
                              ont           = "BP",
                              OrgDb         = "org.Dm.eg.db",
                              pvalueCutoff  = 1,
                              qvalueCutoff  = 0.1,
                              pAdjustMethod = "BH")
dotplot(compPathway, showCategory = 10, title = "GO Enrichment Analysis") + theme(axis.text.y = element_text(size=8))

compPathway <- compareCluster(geneCluster   = genes,
                              fun           = "enrichKEGG",
                              keyType       = "ncbi-geneid",
                              organism      = "dme",
                              pvalueCutoff  = 1,
                              qvalueCutoff  = 0.2,
                              pAdjustMethod = "BH")
dotplot(compPathway, showCategory = 10, title = "KEGG Enrichment Analysis") + theme(axis.text.y = element_text(size=8))

ISC

genes = lapply(list(c1=isc_c1, c2=isc_c2, c3=isc_c3, c4=isc_c4, c5=isc_c5, c6=isc_c6, c7=isc_c7), function(i){ensembl.genes[rownames(i)]$entrezgene_id})
compPathway <- compareCluster(geneCluster   = genes,
                              fun           = "enrichGO",
                              keyType       = "ENTREZID",
                              ont           = "BP",
                              OrgDb         = "org.Dm.eg.db",
                              pvalueCutoff  = 1,
                              qvalueCutoff  = 0.1,
                              pAdjustMethod = "BH")
dotplot(compPathway, showCategory = 8, title = "GO Enrichment Analysis") + theme(axis.text.y = element_text(size=8))

compPathway <- compareCluster(geneCluster   = genes,
                              fun           = "enrichKEGG",
                              keyType       = "ncbi-geneid",
                              organism      = "dme",
                              pvalueCutoff  = 1,
                              qvalueCutoff  = 0.1,
                              pAdjustMethod = "BH")
dotplot(compPathway, showCategory = 10, title = "KEGG Enrichment Analysis") + theme(axis.text.y = element_text(size=8))

Session Info

sessionInfo()
## R version 4.1.2 (2021-11-01)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] grid      stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] ReactomePA_1.38.0                     enrichplot_1.14.1                    
##  [3] clusterProfiler_4.2.0                 png_0.1-7                            
##  [5] ggpubr_0.4.0                          enrichR_3.0                          
##  [7] DT_0.19                               VennDiagram_1.7.0                    
##  [9] futile.logger_1.4.3                   RColorBrewer_1.1-2                   
## [11] org.Dm.eg.db_3.14.0                   scales_1.1.1                         
## [13] reshape2_1.4.4                        knitr_1.36                           
## [15] biomaRt_2.50.0                        GenomicFeatures_1.46.1               
## [17] AnnotationDbi_1.56.2                  genefilter_1.76.0                    
## [19] DESeq2_1.34.0                         SummarizedExperiment_1.24.0          
## [21] Biobase_2.54.0                        MatrixGenerics_1.6.0                 
## [23] matrixStats_0.61.0                    BSgenome.Dmelanogaster.UCSC.dm6_1.4.1
## [25] BSgenome_1.62.0                       rtracklayer_1.54.0                   
## [27] GenomicRanges_1.46.0                  Biostrings_2.62.0                    
## [29] GenomeInfoDb_1.30.0                   XVector_0.34.0                       
## [31] IRanges_2.28.0                        S4Vectors_0.32.2                     
## [33] BiocGenerics_0.40.0                   forcats_0.5.1                        
## [35] stringr_1.4.0                         dplyr_1.0.7                          
## [37] purrr_0.3.4                           readr_2.0.2                          
## [39] tidyr_1.1.4                           tibble_3.1.6                         
## [41] tidyverse_1.3.1                       EnhancedVolcano_1.12.0               
## [43] ggrepel_0.9.1                         ggplot2_3.3.5                        
## [45] pheatmap_1.0.12                       cluster_2.1.2                        
## 
## loaded via a namespace (and not attached):
##   [1] utf8_1.2.2               tidyselect_1.1.1         RSQLite_2.2.8           
##   [4] htmlwidgets_1.5.4        BiocParallel_1.28.0      scatterpie_0.1.7        
##   [7] munsell_0.5.0            withr_2.4.2              colorspace_2.0-2        
##  [10] GOSemSim_2.20.0          filelock_1.0.2           highr_0.9               
##  [13] ggalt_0.4.0              rstudioapi_0.13          ggsignif_0.6.3          
##  [16] DOSE_3.20.0              Rttf2pt1_1.3.9           labeling_0.4.2          
##  [19] GenomeInfoDbData_1.2.7   polyclip_1.10-0          farver_2.1.0            
##  [22] bit64_4.0.5              rprojroot_2.0.2          downloader_0.4          
##  [25] treeio_1.18.0            vctrs_0.3.8              generics_0.1.1          
##  [28] lambda.r_1.2.4           xfun_0.29                BiocFileCache_2.2.0     
##  [31] R6_2.5.1                 graphlayouts_0.7.1       ggbeeswarm_0.6.0        
##  [34] locfit_1.5-9.4           gridGraphics_0.5-1       bitops_1.0-7            
##  [37] cachem_1.0.6             fgsea_1.20.0             DelayedArray_0.20.0     
##  [40] assertthat_0.2.1         BiocIO_1.4.0             ggraph_2.0.5            
##  [43] beeswarm_0.4.0           gtable_0.3.0             ash_1.0-15              
##  [46] tidygraph_1.2.0          rlang_0.4.12             splines_4.1.2           
##  [49] lazyeval_0.2.2           rstatix_0.7.0            extrafontdb_1.0         
##  [52] checkmate_2.0.0          broom_0.7.10             yaml_2.2.1              
##  [55] abind_1.4-5              modelr_0.1.8             crosstalk_1.2.0         
##  [58] backports_1.3.0          qvalue_2.26.0            extrafont_0.17          
##  [61] tools_4.1.2              ggplotify_0.1.0          ellipsis_0.3.2          
##  [64] jquerylib_0.1.4          Rcpp_1.0.7               plyr_1.8.6              
##  [67] progress_1.2.2           zlibbioc_1.40.0          RCurl_1.98-1.5          
##  [70] prettyunits_1.1.1        viridis_0.6.2            cowplot_1.1.1           
##  [73] haven_2.4.3              fs_1.5.0                 magrittr_2.0.1          
##  [76] data.table_1.14.2        futile.options_1.0.1     DO.db_2.9               
##  [79] reactome.db_1.77.0       reprex_2.0.1             patchwork_1.1.1         
##  [82] hms_1.1.1                evaluate_0.14            xtable_1.8-4            
##  [85] XML_3.99-0.8             readxl_1.3.1             gridExtra_2.3           
##  [88] compiler_4.1.2           maps_3.4.0               shadowtext_0.0.9        
##  [91] KernSmooth_2.23-20       crayon_1.4.2             htmltools_0.5.2         
##  [94] ggfun_0.0.4              tzdb_0.2.0               aplot_0.1.1             
##  [97] geneplotter_1.72.0       lubridate_1.8.0          DBI_1.1.1               
## [100] tweenr_1.0.2             formatR_1.11             dbplyr_2.1.1            
## [103] proj4_1.0-10.1           MASS_7.3-54              rappdirs_0.3.3          
## [106] Matrix_1.3-4             car_3.0-12               cli_3.1.0               
## [109] igraph_1.2.8             parallel_4.1.2           pkgconfig_2.0.3         
## [112] GenomicAlignments_1.30.0 xml2_1.3.2               ggtree_3.2.0            
## [115] annotate_1.72.0          vipor_0.4.5              bslib_0.3.1             
## [118] rvest_1.0.2              yulab.utils_0.0.4        digest_0.6.28           
## [121] graph_1.72.0             rmarkdown_2.11           cellranger_1.1.0        
## [124] fastmatch_1.1-3          tidytree_0.3.5           restfulr_0.0.13         
## [127] curl_4.3.2               graphite_1.40.0          Rsamtools_2.10.0        
## [130] rjson_0.2.20             nlme_3.1-153             lifecycle_1.0.1         
## [133] jsonlite_1.7.2           carData_3.0-4            viridisLite_0.4.0       
## [136] fansi_0.5.0              pillar_1.6.4             lattice_0.20-45         
## [139] ggrastr_1.0.1            KEGGREST_1.34.0          fastmap_1.1.0           
## [142] httr_1.4.2               survival_3.2-13          GO.db_3.14.0            
## [145] glue_1.5.0               bit_4.0.4                ggforce_0.3.3           
## [148] stringi_1.7.5            sass_0.4.0               blob_1.2.2              
## [151] memoise_2.0.0            ape_5.5